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ABSTRACT 

Accurate relative spectrophotometry is critical for many science applications. Small wavelength scale 
residuals in the flux calibration can significantly impact the measurements of weak emission and 
absorption features in the spectra. Using Sloan Digital Sky Survey data, we demonstrate that the 
average spectra of carefully selected red-sequence galaxies can be used as a spectroscopic standard to 
improve the relative spectrophotometry precision to 0.1% on small wavelength scales (from a few to 
hundreds of Angstroms). We achieve this precision by comparing stacked spectra across tiny redshift 
intervals. The redshift intervals must be small enough that any systematic stellar population evolution 
is minimized and less than the spectrophotometric uncertainty. This purely empirical technique does 
not require any theoretical knowledge of true galaxy spectra. It can be applied to all large spectroscopic 
galaxy redshift surveys that sample a large number of galaxies in a uniform population. 
Subject headings: galaxies:evolution - methods: data analysis — galaxies: emission lines — galaxies: 
absorption lines — techniques: spectroscopic 
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1. MOTIVATION 

In astronomical spectroscopy, accurate spectrophoto- 
metric calibration is critical for many science applica- 
tions, but difficult to achieve. All calibrations require 
a standard source whose intrinsic properties we know. 
By observing it with the same system (atmosphere, tele- 
scope, instrument, detector) as used for observing the sci- 
ence targets, we can infer the throughput of the system as 
a function of wavelength. However, the accuracy of this 
calibration is limited by how well we know the standards 
and their measurement uncertainty. For ground-based 
astronomical obser vations, we have to rely on natural 
standards (but see iKaiser et all 120081 . 120101 and lAlbertl 
120111 for attempts to use satellite-mounted standards). 
Nearly all extragalactic spectroscopy ha s used stars a s 
the calibration standards (e.g.. Stougli ton et al.|[2002l ). 
Therefore, the accuracy is limited by how well we know 
the true spectra of stars and the uncertainty in the stan- 
dard star spectra. 

The current state-of-the-art spectrophotometric cal- 
ibration utilized in large surveys c an be found in 
the Sloan Digital Sky Survey (SDSS; lYork et al.ll2000l: 
lAbazaiian etal.1 [20091 ). In SDSS, each plate includes a 
set of 16 standard stars, which are color-selected to be F8 
subdwarfs. Their spectra are reduced in the same way 
as science targets. The calibration is achieved by com- 
paring the galactic-reddening-corrected spectra of these 
stars to a grid of theoretical spectra generated from Ku- 
rucz model atmospheres. The average ratio of each star 
to its best-fit model is taken as the flux c alibration vector 
and applied to all the science targets (jAbazaiian et al.l 
|2004[ ). Finally, the corrected spectra are put on an ab- 
solute scale by comparing the synthetic photometry to 
fiber magnitudes from the SDSS photometry. 

The accuracy of this method is limited by the signal 
to noise of the standard star spectra and the accuracy 
of the models. The former limits the confidence in us- 
ing them to correct for the system response variation on 
small scales. In SDSS, a fourth-order b-spline with 50 



breakpoints is used to smooth the flux calibration vec- 
tor before applying it to calibrate science targets, thus 
leaving small-scale variations uncorrected. The theoreti- 
cal model of the standard stars could also have unknown 
systematic residuals on small scales. 

Overall, the wavelength-dependent relative flux cali- 
bration in SDSS is good to the level of l%-2%. However, 
this accuracy is still inadequate in certain science appli- 
cations. The problem is illustrated in Figure [T] Here, we 
show the [N II] /[S II] flux ratio and i:»„(4000) index for a 
population of galaxies that have emission lines with ra- 
tios characteristic of low- ionization nuclear emission-line 
regions (LINERs). The sample is volume-limited. Their 
[NII]/[SII] ratios vary systematically with redshift, at 
a level nearly comparable to the intrinsic scatter at a 
single redshift. Similarly, the D„(4000) index also show 
unphysical variations with redshift. 

These emission line fluxes and spectral indices are mea- 
sured with our own software, which is an impr oved ver- 
sion of the code described bv lYan et al.l (|2006[ ). Similar 
systematics also exist if we adopt the D R7 catalog pro- 
duced by the MPA-JHU collaboratiorQ (i Tremonti et all 
l2004f ) , as shown in the right panels of Figure [TJ The de- 
tailed wiggles differ between our measurements and the 
MPA/JHU catalog. These differences are mainly due to 
the difference in the methods used to measure the con- 
tinuumo However, the root cause of the systematic vari- 
ation with redshift in both cases is due to a small-scale 
(~ lOA) systematic flux calibration residual, as we will 
demonstrate in this paper. 

In this paper, we will establish a new method for im- 
proving the flux calibration accuracy by an order of mag- 
nitude using coadded red galaxy spectra. Because red 
galaxies have very uniform spectra and evolve very slowly 
with redshift, when coadded they provide a high signal- 

^ http://www.mpa-garching.mpg.de/SDSS/ 

^ Our code measures the continuum in the sidebands bracketing 
the emission line, while the MPA/JHU code measures it using a 200 
pixel median smoothing of the emission-line-subtracted continuum. 
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Fig. 1.— Left panels show the distribution of the [N II] A6583/[S II] AA6716,6731 flux ratio (upper panel) and D„(4000) (lower panel) 
for LINERs as a function of redshift. The sample is volume-limited with Mo.i^ < —20.7. LINERs are selected using the criteria described 
bv lKewlev et all ^2006), except that no criterion involving [O I] A6300 is required. The dark points with tiny error bars show the median 
of the distribution. The error bar indicates the error of the median, computed using the formular given by Beers et al. ( 1990). The wiggly 
systematic variations with redshift are unphysical. The right panels are similar to the left but using measurements from the MPA/JHU 
DR7 catalog. The exact measurement methods differ but they all show significant systematic variations with redshift, which are due to a 
small-scale systematic flux calibration residual. 



to-noise standard. In addition, because this standard can 
be found at any redshift, no wavelength will be at disad- 
vantage in the calibration due to features in the standard. 
We will show that this method can achieve a relative flux 
calibration accuracy of 0.1% on scales of a few hundred 
of Angstroms in SDSS. 

In Section 2, we show how to measure the flux cali- 
bration residual in the SDSS using stacked red galaxy 
spectra, and test the improved spectrophotometry. We 
discuss potential improvements in Section 3 and summa- 
rize in Section 4. 

2. IJVIPROVED SPECTROPHOTOIMETRY 

We employ the SDSS Data Release 7 IMain galaxy sam- 
ple from the New York Univ ersity Value-Added Galaxy 
Catalog (jBlanton et al.|[2"005D . The absol ute m agnitudes 
are derived using Blanton & Roweis's (|2007D k-correct 
codes (version 4_2). 

The first step in our calibration procedure is to identify 
red galaxies with old stellar populations. We first apply 



two cuts 
galaxies: 
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color to select the red sequence 



°-\g-r) > -Om{Mo.i., 
"■\g-r) < -0.02(Mo.v 



51og/i) + 0.49 (1) 

5 log /i) -I- 0.59, (2) 

where all magnitudes are in the AB system. The su- 
perscript 0.1 indicates that the filter systems are shifted 
blueward by a factor of 1.1. 

We bin all red galaxies with 0.06 < z < 0.15 that 
satisfy the above criteria and have apparent r-band 
model magnitude brighter than 17.7 into 118 narrow red- 
shift bins. The binsize corresponds to 3 pixels in the 
SDSS wavelength grid, which is logarithmically spaced 



(A log A = 10 "*). This spacing translates to Alog(l + 
z) = 3 X 10"''. We do not need to apply any absolute 
magnitude limit to our sample, because we will only com- 
pare each redshift bin with its neighbors, and the redshift 
differences between adjacent bins are tiny (~ 0.0008). 
Essentially, we are always comparing samples with the 
same luminosity distribution. 

Some red-sequence galaxies have star formation but 
appear red due to dust extinction. To make our standard 
as stable as possible, we exclude, in each redshift bin, 
the 20% of red sequence galaxies that have the lowest 
Z3„(4000) measurement. The final number of galaxies in 
each redshift bin ranges from 700 to 1300. 

We first de-redden the spec tra for galactic extinction 
according t o the dust map of iSchlegel et al.l (J1998.) and 
lO'Donnelil ()1994[ ) extinction curve. Each spectrum is 
then de-redshifted and smoothed to the same resolution 
of 350kms" by convolving it with a Gaussian kernel with 
varying width, determined according to the fixed instru- 
mental resolution and the measured velocity dispersion 
of each galaxy. We normalize each spectrum so that they 
all have the same median flux between 6010A and 6100A. 
Finally, we average them in each redshift bin to obtain 
the stacked spectrum. 

The top panel of Figure [5] shows two stacked spectra 
from a pair of adjacent redshift bins near z '^ 0.1. There 
are very tiny, nearly invisible differences between these 
two coadds. The ratio between the two stacks (bottom 
panel) varies with wavelength slightly, deviating from 
1 by a few percent. The differences have two possible 
sources. One possibility is intrinsic difference between 
the two populations, which is a function of rest-frame 
wavelength. The other possibility is the fiux calibration 
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Fig. 2. — Top panel shows a comparison of two stacked spectra of rcd-scquence galaxies from adjacent narrow redshift bins near z ~ 0.1 
(one solid black curve and one dashed red curve). The two spectra overlap to very high precision. Bottom panel shows the ratio between 
the two spectra. The tiny differences between the stacked spectra are dominated by spectrophotometric calibration residuals, as the same 
rest-frame features are redshifted to slightly different wavelengths. The intrinsic differences between the two populations are subdominant 
as demonstrated below. 



residual, since the same spectral features are redshifted 
to slightly different wavelengths. The latter is a func- 
tion of the observed wavelength. We can check which 
source dominates the difference by comparing the ratios 
from multiple pairs of stacks in the rest frame and in the 
observed frame. 

In Figure |3l we take the ratios from 117 pairs of 
stacked spectra in adjacent redshift bins, normalize them 
at 6100A, and plot them against the rest-frame wave- 
length. (We always divide the higher-z bin stack by the 
lower-z bin stack in each pair.) Lining them up in rest 
frame smears out the flux calibration difference but high- 
lights the intrinsic population difference. There is very 
little consistent variation with wavelength among these 
ratio curves. The median ratio has tiny deviations from 
1. Over the wavelengths between 3553. SA and 8072. 4A 
which are completely covered by all 117 ratio curves, the 
standard deviation over all pixels in the median curve 
is 0.00027; only 1.3% of all pixels have a median ratio 
deviating from 1 by more than 0.1%. 

Although the median curve does not show any sig- 
nificant systematic variation, each ratio curve still con- 
tains the intrinsic difference between the two redshift 
bins, which could be due to a tiny mismatch in galaxy 
properties between the two populations. This can be 
seen by the larger variation at particular wavelengths, 
such as [O II] A3727, Mg b A5170, Na D A5890,5896, Ha, 
[N II] A6583, etc. The intrinsic difference could also give 
rise to a small smooth deviation on large scales. For ex- 
ample, the ratio curve in the bottom panel of Figure [2] 
deviates systematically from 1 on the blue side and has a 
large-scale downward slope. Because these intrinsic dif- 
ferences happen in both positive and negative directions 
with roughly equal probabilities, they largely cancel out 
in the median curve. We could also remove the large- 
scale component by smoothing in the following analysis. 

The generally larger variation on the blue end com- 
pared to the red end could also be due to the lower flux 
and thus smaller signal-to-noise ratios in the blue. This 



also partly explains the larger variation seen around Ca II 
H and K lines. 

In Figure m we plot the same ratio curves against ob- 
served wavelength. This smears out the intrinsic popu- 
lation difference but highlights the flux calibration resid- 
uals. Here, we divide each ratio curve by its 500 pixel 
median-smoothed version to remove any large scale com- 
ponent. Any small-scale features due to intrinsic popula- 
tion difference will be smeared out by redshifts and can- 
celled out among the multiple pairs. Thus, they do not 
introduce any systematic variation as a function of the 
observed wavelength. Evidently, the contribution from 
the flux calibration residual is significant at the level of 
1%, much larger than the median deviation due to pop- 
ulation differences (< 0.1%). Since our redshift bins are 
3 pixels wide, the ratio between two stacked spectra at 
wavelength Afc (fc is the pixel number) is the ratio of the 
system response function (i?(A)) between the two 3 pixel 
wide elements centered on Afe+3 and Xk, respectively. 

Figure|3]and|4]demonstrate that the difference between 
these stacked spectra can be separated into two main 
components: a systematic component due to flux cali- 
bration error and a noise-like component due to small 
intrinsic population differences. The latter is only signif- 
icant blueward of 5000A. With many pairs of adjacent 
redshift bins to average out the "random" population 
differences, the stacked spectra of old red galaxies across 
small redshift intervals can be used as a spectroscopic 
standard to constrain relative system response variations 
on small scales. 

To express the above procedure in mathematical terms, 
we use /(A) to denote the average spectrum of red galax- 
ies which has been shown to display only tiny systematic 
variation with redshift. The observed stacked spectrum 
at redshift z is denoted as g{X). We use i?(A') to de- 
note the to-be-constrained system response function in 
the observed wavelength space (A' = A(l + zj). Since 
we coadded spectra in 3 pixel wide redshift bins, we are 
effectively using a 3 pixel box smoothed version of R, 
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Fig. 3. — Distribution of ratios between pairs of stacked spectra from adjacent redshift bins as a function of rest-frame wavelength. There 
are 117 ratio curves from 118 stacked spectra. The thick black line shows their median, and the thinner lines represent the 5, 25, 75, and 
95 percentiles at each pixel. The systematic difference among these stacked spectra is tiny, which demonstrates that they can be used as 
spectroscopic standards to calibrate spectrophotometry. See Figure |4] 
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Fig. 4. — Distribution of ratios between pairs of stacked spectra from adjacent redshift bins as a function of observed wavelength. There 
are 117 ratio curves from 118 stacked spectra. The thick black line shows their median, while the thinner lines represent the 5, 25, 75, and 
95 percentiles at each pixel. Clearly, there is systematic variation with observed wavelength. Comparing this with Figure |3] shows that the 
spectrophotometry calibration error dominates the differences between the stacked speci;ra. Note the tiny scale on the vertical axis. This 
shows we can achieve high- precision spectrophotometry with stacked spectra of red-sequence galaxies. 



which we denote as R{X'). 
For two adjacent redshift bins, Zi and z^+i, we have 

ft(A) = i?(A(l + z,))/,(A) (3) 

g,+i(A) = i?(A(H-z,+i))/,+i(A) (4) 

We define the ratio, Qi^i+i — gij^i/gi. Assuming fi = 



fi+i, we have 

hi.R(A(l-f z,+i)) 



lni?(A(l + z,))=ln9M+i(A) (5) 



To rewrite this equation in the observed wavelength 
space, we set 

Afe(l + z,) = A'fe, (6) 

where k is the pixel index, fc = 0, l,...,n — 1, and n is 
the total number of pixels. Since the wavelength grid is 
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logarithmically spaced, and the spacing between adjacent 
redshift bins is equivalent to 3 pixels, 

Afc(l + zj+i) = Afc+3. (7) 

We shift each qi,i+i into the observed frame and define 

qU+liK)^<l^,^+l{>^k)■ (8) 

Substituting Equations dUl-dH) into Equation ([5]), we 
have _ _ 

lni?(A;,+3) - lni?(A;,) = hiql^^.iX',), (9) 

where fc = 0, 1, . . . , n — 4; i = 0, 1, . . . , rn — 2; and m is 
the total number of redshift bins. 

Subtracting the smooth component in each q'^ ,j_|_]^ and 
taking their median at each A'j,, we have, 

lni?(A',+3) -lni?(A;) = Ing-'(A^) (10) 

which we abbreviate as 

\nRk+3-\nRk^\nq',,. (11) 

Given this equation array, we would like to solve for R 
using the measured q' . With n unknowns but only n — 3 
equations, we are short of three boundary conditions. In- 
accurate boundary conditions will introduce oscillations 
in the solution. Here, we define Ini? to be the 3 pixel 
box smoothed version of Ini?, i.e.. 



\nRk = (lni?fe_i 
It is easy to show that 



lni?fc + lni?fc+i)/3. (12) 



Ini? 



fe+2 



-lni?fe+i = 



In i?fe+3 - In Rk 



= \nq'iK) (13) 



Now it is trivial to solve for In R by setting an ar- 
bitrary boundary condition. It is arbitrary because we 
have no constraint on the absolute_scale of R. We then 
take the exponential and normalize R so that the average 
between 3800 A and 9200 A is 1. Figure [5] shows the final 
spectrophotometry calibration. The correction factor is 
tabulated in Table [TJ To correct spectra in SDSS, divide 
the spectra by the vector listed here. 

The median q' has a median uncertainty of 0.00009 
over all pixels, which is the median uncertainty on the 
response difference between adjacent 3-pixel-wide ele- 
ments. Therefore, the final R will have an accumulated 
error of ~ 0.0009 over every 300 pixels. Therefore, we 
have achieved 0.1% or better accuracy in relative spec- 
trophotometry on wavelength scales shorter than 300- 
600A. 

This uncertainty is of the same order as the uncertainty 
propagated from the error in the stacked spectra. The 
latter is limited by the number of galaxies in the stack 
and the signal-to-noise ratio of each spectrum. Increas- 
ing either the signal-to-noise of individual spectra, the 
number of galaxies in each bin, or the total number of 
redshift bins can improve the precision of this calibra- 
tion. The uncertainty is also larger on top of sky lines 
(atmosphere emission lines, e.g., 5577A, 6300A) than in- 
between sky lines due to the lower signal-to-noise ratio 
at those wavelengths and the non-Gaussian nature of the 
sky line subtraction residuals. The larger uncertainties 



at short wavelengths (see Figure S]) are due to the lower 
signal-to-noise ratios of red galaxy spectra toward the 
blue end. 

The original flux calibration in SDSS is already ac- 
curate to l%-2% level. Now we can remove these few 
percent residuals to make it accurate to subpercent level 
by using red-sequence galaxies as the spectroscopic stan- 
dard. 



TABLE 1_ 

Residual System Response (R) as a Function of 

Wavelength 



Vacuum Wavelength 


(A) 


R 


3776.5914 




1.01105 


3777.4611 




1.01877 


3778.3310 




1.02491 


3779.2011 




1.02918 


9238.4689 




1.00488 


9240.5964 




1.00352 


9242.7243 




1.00060 


9244.8528 




0.99609 


9246.9817 




0.99077 



(This table is available in its entirety in machine-readable and Vir- 
tual Observatory (VO) forms in the online journal. A portion is 
shown here for guidance regarding its form and content.) 

In this residual correction curve, one can see absorp- 
tion features at 3933A, 3968A, 5890A, and 5896A, which 
could be due to Ga II and Na I absorption in the Milky 
Way's interstellar medium (ISM), or maybe the slightly 
incorrect F star models around these absorption features. 
There are also features which line up with small scale 
features in the atmospheric transparency curve, such as 
the water absorption bands around 7150-7400A, 8100- 
8400A, and above 8900A. 

To test our improved flux calibration, we applied the 
correction derived above to the SDSS spectra and remea- 
sured the emission lines for all galaxies in DR7. In the 
right panel of Figure [HI we show the [N II] /[S II] ratios 
and i?„(4000) measurements for LINERs selected in the 
same manner as in Figure [TJ The latter is reproduced 
in the left panel here. These are largely the same galax- 
ies, except that ^ 7% of them differ du e to ch anges in 
their classifications in the iKewlev et al.l (|2006[ ) scheme. 
This indicates that the impact of the flux calibration on 
the classification is not negligible. Gomparing the two 
panels, the unphysical systematic trend with redshift 
is largely gone with the new calibration. The median 
[NII]/[SII] ratio becomes fairly smooth with redshift. 
The smooth variation with redshift could result from the 
aperture effect since the same fiber size corresponds to 
different physical scale at different redshifts. The inves- 
tigation of these apert ure-dependent line ratio variations 
(|Yan fc BlantonI [20110 is precisely what motivated us to 
improve the small-scale calibration. Though our calibra- 
tion is still not perfect as shown by the small wiggles 
left in the right panel of Figure [51 it is sufficient for our 
science project. 

3. POTENTIAL IMPROVEMENTS 

Further improvements of the above calibration are pos- 
sible along several paths. 
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Fig. 5. — Residual system response as a function of observed wavelength. The normalization is arbitrary. 
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Fig. 6.— Distribution of [N II] A6583/[S II] AA6716,6731 flux ratio (upper panel) and D„(4000) (lower panel) for LINERs as a function of 
redshift before (left panel) and after (right panel) applying our flux calibration. The unphysical variations with redshift largely disappear 
after applying our flux calibration. (The apparent underdensity of points at ^ ~ 0.114 is caused by the fact that, at this redshift, the large 
sky subtraction residual at 5577A leads to uncertain measurements of [O III] A5007, whose significant detection is required for the selection 
of LINERs.) 



First, we could iterate the above calibration process. 
Since we are measuring a smoothed version of i?, the 
original R must have sharper features than what we mea- 
sure. Thus, folding in a first-order correction and iter- 
ating would help us measure those sharp features more 
accurately. 

Second, we could use finer redshift bins to get the pixel- 
to-pixel variation. However, the advantage of higher 
pixel resolution needs to be weighted against the poorer 
signal-to-noise in stacked spectra from the smaller num- 
ber of galaxies in each bin. 



Third, as shown in Figure [3l certain rest-frame wave- 
lengths in the spectra show larger variation than others, 
such as [O II] A3727, Ca II H and K, Na I D, Ha, etc. 
One way to improve our calibration is to exclude these 
regions (or downweight them) when taking the median 
(or average) in the observed wavelength space. 

Fourth, the calibration residual has features that might 
originate from three sources: the atmosphere, the instru- 
ment, and the ISM. These features could have intrinsic 
variations with different dependencies. For the atmo- 
spheric features, by grouping galaxies according to the 
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observing conditions, such as airmass and/or humidity, 
and doing the cahbration separately, we could reduce the 
intrinsic scatter at these wavelengths and improve the 
calibration. For features that originate from the tele- 
scope and the instrument, it could differ between the 
two spectrographs used in SDSS. It could also vary with 
the ambient environment during observation or with time 
due to the degrading of mirror coatings. Again, group- 
ing galaxies according to the spectrograph, the dates of 
observation, or ambient environment may identify these 
features and achieve better corrections. For the absorp- 
tion in the ISM, galaxies along different lines of sight 
should have different corrections. By grouping red galax- 
ies according to their Galactic extinction, we can get 
better control on these wavelengths and may learn some- 
thing about the ISM at the same time. Certain stud- 
ies may prefer to have the interstellar absorption un- 
corrected. However, we first need to identify robustly 
the source of the absorption. The experiment proposed 
above would also be useful in this respect. 

Finally, the procedure we have described above does 
not constrain at all the smooth, large-scale variation in 
the flux calibration. To constrain that, we could adopt 
wider redshift bins, a larger step size, and better con- 
trol for a uniform population of galaxies at all redshifts. 
In addition, we also need to take into consideration the 
evolution in the average spectra. It is still possible, since 
the evolution is a function of the rest-frame wavelength 
and the spectrophotometry is a function of the observed 
wavelength. Thus, they can be separated. With a large 
enough sample, it may be possible to solve for evolution 
and spectrophotometry simultaneously. 

4. SUMMARY 

We have demonstrated that the stacked spectrum of 
carefully selected red-sequence galaxies can be used as a 
spectrophotometric standard to produce an accurate rel- 
ative flux calibration for large spectroscopic surveys. For 
SDSS, we have achieved the accuracy of ~ 0.1%. This 
method can be applied to surveys like the Baryon Os- 
cillat ion Spectroscopic Survey (]30SS; lEisenstein et al.l 
1201 It ), as long as the survey covers a sufficiently large 
number of red galaxies. In principle, other types of galax- 
ies could also be utilized provided stacking a large enough 



sample can produce a reasonably stable spectrum over 
some wavelengths, which can be verified empirically, as 
done in Figure |3l 

However, this method is not meant to replace the tra- 
ditional calibration methods using standard stars. To 
obtain best calibration, we still recommend the practice 
of observing standard stars simultaneously with science 
targets, as done in SDSS. They are essential for getting 
the smooth component in the spectrophotometry and for 
correcting for those extinction sources that vary on short 
timescales, such as the molecular absorption in the at- 
mosphere. Our method is meant to be complementary 
and to improve a reasonably well calibrated spectropho- 
tometry to the accuracy required for certain science ap- 
plications, such as th e study of weak lin e emission in red 
galaxies, Lya forest (|Slosar et al.ll201ll ). and any studies 
requiring high accuracy spectral index measurements. 
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